Unquenched complex Dirac spectra at nonzero chemical potential: two-colour QCD 

lattice data versus matrix model 
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We compare analytic predictions of non-Hermitian chiral random matrix theory with the complex 
Dirac operator eigenvalue spectrum of two-colour lattice gauge theory with dynamical fermions at 
nonzero chemical potential. The Dirac eigenvalues come in complex conjugate pairs, making the 
action of this theory real, and positive for our choice of two staggered flavours. This enables us to 
use standard Monte-Carlo in testing the influence of chemical potential and quark mass on complex 
eigenvalues close to the origin. We find an excellent agreement between the analytic predictions 
and our data for two different volumes over a range of chemical potentials below the chiral phase 
transition. In particular we detect the effect of unquenching when going to very small quark masses. 

PACS numbers: 12.38.Gc, 02.10.Yn 



Non-Hermitian operators appear in many different ar- 
eas of Physics: in S-matrix theory with absorption 
in neural network and dissipative quantum dynamics 
3, disordered systems with imaginary vector potential 
3j, Quantum Chromodynamics (QCD) with a non-zero 
chemical potential /i 0| , or 0- vacuum term Q . Studying 
the influence of these non-Hermiticity parameters may 
serve to a better understanding of fundamental problems 
in nature, such as the deconfinement transition observed 
in heavy-ion collisions, or the strong CP problem. In 
the present work we use an analytically solvable random 
matrix model (MM) as an effective theory based on sym- 
metries, and focus on the applications to QCD at fi ^ 0. 
Given the wide range of complex MM applications 1] we 
expect our approach to be relevant to other problems in 
the same symmetry class. 

Lattice Gauge Theory (LGT) is one of the most promi- 
nent methods to study nonperturbative QCD. However, 
the introduction of a chemical potential poses a difficult 
problem: the Dirac operator and thus the action becomes 
complex non-Hermitian. In general this invalidates stan- 
dard Monte-Carlo techniques using importance sampling, 
the so-called sign problem. The quenched approximation, 
eliminating the complex fermion determinant, has been 
explained to fail using a MM approach Q. It maps QCD 
to a theory with conjugate quarks, with chiral symme- 
try unbroken in the massless limit (see ^ for a recent 
discussion, including phase quenched MM). 

Since the role of dynamical fermions is crucial our goal 
is to show that MM can predict unquenched lattice data 
at /i 7^ 0. Although analytical MM results exist for un- 
quenched QCD (3, a comparison to LGT remains an 
open challenge, due to a cornplex valued density with 
highly oscillatory regions 0, The various ways to 
avoid the sign problem including reweighting, Taylor ex- 
pansion, imaginary chemical potential or combinations 
of these reviewed in Q rely essentially on being close to 
the chiral phase transition at high temperature T k, Tc. 



These techniques do not apply easily at the region of our 
interest at small T « and where MM are conjectured 
to coincide with a particular limit of the underlying ef- 
fective chiral Lagrangian, the e-regime (exPT). MM are 
a different realisation of the group integral over constant 
Goldstone modes at /i 7^ for QCD |^ ^-nd they 
clearly distinguish the three different chiral symmetry 
breaking (xSB) patterns at ^ 7^ 0. Thus they may 
serve as a test for algorithms solving the sign problem. 
In order to access in this region with dynamical fermions 
we follow the suggestion [lj|: gauge group SU{2) or the 
adjoint representation lead to a real action as complex 
Dirac eigenvalues come in conjugated pairs. We have 
chosen SU{2) with 2 staggered fermions because in this 
symmetry class MM predictions are available The 
difhculty is that large masses m effectively quench the 
Dirac eigenvalues close to the origin. Unquenching is 
seen only for small m depending on the particular size of 
the system. Our preliminary results for SU{2) have been 
presented in compared to SU{2) at ^ = 0| and 
quenched QCD at /i 7^ [3, E3 • 

We start by introducing the MM j the complex ex- 
tension of the chiral Gaussian Symplectic Ensemble, and 
its predictions. For an arbitrary number of Nf quark 
flavours with masses rrif it is defined as: 



Zy{{mf};^i) = /d<I>d*]^det 

U(« + ,.)XN / = 1 

X exp [-A^Tr($t$ + 5»t^)] _ (1) 
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The matrices $ and \1/ of size {N + v) x N contain quater- 
nion real matrix elements without being self-dual, 1 is 
the quaternion unity element. The Gaussian average of 
variance a= I/V2 replaces the gauge average, where we 
have fixed topology or equivalently the number of ex- 
act zero eigenvalues v > 0. The total number of eigen- 
values 2N ^ V relates to the volume. For = 1 this 
MM was shown ll] to be in the symmetry class of the 
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adjoint representation, or gauge group SU(2) with stag- 
gered fermions, as in our case. In the qualitative 
differences between the 3 xSB classes were compared nu- 
merically. Assuming universality, which is true for the 
MM of QCD HHEa, we choose the ^i-part ^ 1, with 
the same symmetry as the kinetic part. This choice in 
eq. permits an eigenvalue representation with Jaco- 



and mass term Jl/'ii N 



(here z = 



iy is real 



at fi — 0). Due to the Jacobian complex eigenvalues are 
repelled both from the x- and y-axis, the signature of this 
xSB class. The QCD MM 7] has no such repulsion. 

In 13] all complex eigenvalue correlations were derived 
in a closed form as the Pfafhan of a matrix kernel of skew 
orthogonal polynomials on C. We only give the density 
for Nf = 2 of equal mass needed here. For a discussion of 
matching MM and staggered flavours we refer to the first 
ofRef. JJJ. No 2n(i or 4t/i root is taken here. Taking the 
large- A'^ limit at weak non-Hermiticity 1] the eigenvalues, 
masses and fi are rescaled as 



z-2Na EE ^, m-2Na = ri, 



N- 



hm fi^-2N : 



(2) 



resulting in the usual MM level spacing tt. The e-regime 
of the corresponding chiral Lagrangian dictates the 
same scaling with two different constants: VTim and 
F^V^^, the chiral condensate S and decay constant F-^, 
respectively. This implies the existence of one free fit pa- 
rameter in in our model. This observation has very 
recently been proposed to measure for QCD using 
isospin chemical potential [l^. We note that MM pre- 
dictions at /x = are parameter-free in units of S. 

In terms of the rescaled variables eq. Q the spectral 
density defined as p{z) = (n7=i S'^^Hz — Zj)) is given by 



p(^/=o)(e) = (r2-e'Me,r)'*(e,r) 

for quenched. The rescaled weight and kernel read 



1 



(3) 



(4) 



X (J2.(2V^0^2.(2\/iC*) O), (5) 



containing K- and J-Bessel functions. The weight for 
the eigenvalues is non-Gaussian due to the non-trivial 
decoupling of eigenvectors. The unquenched density 



K{ir], irj) 



(6) 



with a mass dependent correction term is shown in Fig. ^ 
left (for more details, see Compared to = in 

(isf the peaks locating the eigenvalues are split in two. 
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FIG. 1: The spectral density with iV/ = 2 eq. ® at a = 0.012, 
ry = 8.74 and v = (left), and a cut parallel to the a;-axis for 
different masses (right). 

the particularity of this xSB class. A cut parallel to 
the X-axis along the maxima illustrates the effect of dy- 
namical flavours at increasing 77. At value 77 w 100 the 
smallest eigenvalues are effectively quenched, following 
eq. Increasing a rapidly damps the oscillations. In 
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FIG. 2: The density eq. ® at a = 3.28, = 19, and u = 
(left), and a cut parallel to the j/-axis (right) vs. quenched. 

Fig. 12 the density is now constant along a stripe, as ex- 
pected from mean field. The effect of the additional zeros 
of p^^^^^HC) at ^ = ±177 in the left corners of Fig. [3 are 
clearly visible. Here we choose a cut in j/-direction on 
the maxima to highlight the difference to quenched. At 
this value of a the perpendicular cut in cc-direction shows 
no dependence on 77, in contrast to the previous figure. 
In the strong non-Hermiticity limit a ^ 00 the plateau 
extends to 00 in both x- and 7/-direction • Away from 
the edge Fig. |21is already very close to a = 00, to which 
we compared our preliminary data [l^ . 

We now turn to the LGT side, describing the details 
of our simulations summarised in Tabled Our data were 
generated for gauge group SU{2) with Nf = 2 unim- 
proved staggered fermion flavours of equal mass at cou- 
pling constant f3 = — 1.3, using the code of [l^ . 
In this setup the fermion determinant remains real and 
standard Monte Carlo applies. Our choice of parameters 
was dictated by several constraints: i) Wilson fermions 
explicitly break chiral symmetry and are already complex 
at /i = 0, making the effect of /i ^ difficult to disentan- 
gle. At /i = Ginsparg- Wilson type fermions have been 
shown to exactly preserve chiral symmetry and topology 
(see (23 for references and a comparison to MM). Simi- 
lar results have been obtained with improved staggered 
fermions for QCD [2]| . Apart from being very expensive 
it is conceptually not clear so far how to extend these 
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TABLE I: Summary of simulation parameters at /3 = 1.3 



V 




mass ma 


level spacing d 


n 




i^2/E 


Eth 


no. of confi 


6* 


1.0 ■ IQ-^ 


0.035 


8.0(2) • 10"^ 


13.7 


0.0124(3) 


0.39(3) 


~ 17 


45 000 


6* 


1.0 • IQ-^ 


0.07 


6.6(2) • 10"^ 


33.3 


0.0127(3) 


0.34(3) 


~ 25 


15 000 


6* 


1.0 ■ IQ-^ 


20 


5.0(2) ■ 10"^ 


12 566 


0.0153(3) 


0.37(3) 


> 100 


10 000 


6* 


0.2 


0.06 


9.9(1) ■ 10"^ 


19.0 


2.21(4) 


0.38(1) 


~ 10 


20 000 


6* 


0.2 


20 


4.5(1) ■ 10"^ 


13 962 


3.28(4) 


0.38(2) 


> 200 


20 000 


8* 


5.625 • 10-'' 


0.035 


3.0(3) • 10"^ 


35.5 


0.0107(9) 


0.35(5) 


~ 20 


8 000 


8* 


5.625 • 10"'' 


0.07 


2.1(1) ■ 10"^ 


104.7 


0.0124(3) 


0.32(4) 


~ 35 


5 000 


8* 


5.625 ■ 10"'' 


20 


1.6(1) ■ 10"^ 


39 270 


0.0152(3) 


0.37(4) 


> 35 


3 000 


8* 


0.1125 


20 


1.4(1) ■ 10"^ 


44 880 


3.29(4) 


0.38(5) 


> 40 


6 000 



results to ^ 7^ 0. Thus we use the Dirac operator of 
standard, unimproved staggered fermions: 

^a=.y = ^[ {U^{x)K^{x)6y, -h.c.) (7) S 

j=x,y,z 

+(t/,-(a:)«;,-(a:)e''<5^,,^"~ t/|(2/)AcKy)e-''<5^,,_i)] , 

with Hnk variables U and staggered phases k. Conse- 
quently our simulations are topology blind, and we set 

= in the following. We note that even at = 
the transition from staggered to the correct continuum ^ 
symmetry class has not yet been observed, ii) A proper 
resolution of complex eigenvalues requires of the order 
of 10^ configurations for each parameter set. Even with 
cheap staggered fermions we are confined to small lat- 
tices V = 6'^ and 8*. iii) We want to make the window 
where MM apply as large as possible. At the Thouless 
energy scale |23|, Eth ~ -P't/^V^j the Dirac eigenvalues 
start to feel the propagating Goldstone modes, and leave 
the e-regime. Increasing (3 (decreasing Ft^) at fixed V 
fewer Dirac eigenvalues follow MM statistics. This leads 
us to rather strong coupling, iv) Small masses m consid- 
erably slow down the computation. Because of rescaling 
ri = mT,V this effects our choice of the volume too. 

To compare our data with the prediction eq. © we 
have to determine two parameters from the data. First, 
we measure the mean level-spacing d ~ l/p{0), by us- 
ing the Banks-Casher relation at fi — 0, 7rp(0) = Sy, 
where p(0) is the mean spectral density. We project the 
complex eigenvalues onto the real axis and measure their 
distance. Due to /i ^ 1 the 2d geometric distance be- 
tween eigenvalues agrees within errors. This provides the 
volume rescaling of lattice eigenvalues, masses and /i^: 

za-n/d = S, , ma-ir/d =7], and ii^C-n/d=a^. (8) 

Second, the remaining constant C ^ F^/Y, is obtained by 
a fit to the data. We cut the data rescaled as in eq. © ^ 
parallel to the y-axis on the first maxima, which are easy 
to identify for all a. Then we fit to the integral of the 
analytic curve, thus eliminating the choice of bins in this 
direction (for clarity we show histograms in all Figures). 
This procedure avoids cumbersome 3d fits and fixes the 
normalisation of the first maxima. 
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FIG. 3: The effect of dynamical fermions at = 10""^, 
V = 6*: Q = 0.0127, 77 = 33.3 (upper), and a = 0.0124, 77=13.7 
(lower) . 

We start with a <C 1. For mass 77 = 33.3 in Fig. El 
(up) the data show an excellent match to the predic- 
tion eq. (jnj including the 4t/i maximum, clearly deviat- 
ing from the quenched curve shown for comparison. The 
perpendicular cut on the first maxima (up right) is per- 
fectly agreeing, allowing to fit a to an relative error of 
w 2%. The repulsion from the a;— axis correctly identifies 
this xSB class (QCD has a single peak, see and we 

have observed this pattern down to ^= 10"^. Note that 
fitting a does not move the maxima in the left curves. 
A similar data set for y = S** leads to the same conclu- 
sions (see Table 1). Fig. El (down) shows the smallest 
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FIG. 4: The scaling ^i^Vi =const.: l/i=6'', = 10-^ (green) 
vs. V2=8'', /i2 = 5,625-10-* (red), both at ma = 20 for y > 0. 
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FIG. 5: The effect of dynamical fermions at /i = 0.2, V = 6*: 
a = 2.21, r?=19.0 for y > 0. Note that Eth^W. 

mass we could simulate, with 77= 13.7. The data follow 
the analytic curve including the 2nd minimum, deviating 
substantially from quenched. The onset of Exh ~ 17 is 
responsible for a smaller window of agreement (Exh ^ 25 
at 77 = 33.3). The perpendicular cut shows again perfect 
agreement. Despite excessive statistics the fluctuations 
on the maxima are still visible. To exclude that the ob- 
served deviations from quenched are finite volume effects 
we compare to [Tsf at /i = and otherwise identical pa- 
rameters. There finite volume shifts the data to the left 
of the MM curve. Thus for larger volumes the deviation 
from quenched would be even more pronounced. Second 
we check the scaling hypothesis from MM and exPT, 
choosing n\/ ix\ = 6^/8^ with ^fVi = const. For simplic- 
ity [2^ we have chosen a large rj to effectively quench 
the small eigenvalues. In Fig. 01 both sets of data show 
an excellent agreement to the same curve from eq. 
Independent fits for a agree within errorbars. Here only 
?/ > is shown for a better resolution, the second peak 
at ?/ < follows from symmetry. The same results are 
obtained for smaller, unquenched masses (see Table 1). 

Next wc turn to a larger ^, remaining below the phase 
transition ~ 0.3 found in Q for F = 6^, (3 — 1.3. 
Above the transition at /x = 0.4 the data deviate from 
the MM as expected. We first determine d, which fixes 
the scale and rj. At these parameter values unquenching 
is seen in y-direction (Fig. |S1 right), which is used to fit 
a. Here both a and the given 77 influence the width. We 
therefore compare the unquenched curve to the quenched 
one with a' ^ a fitted independently (lower curve). For 
comparison we also give the quenched curve at the same 
a (upper curve). Again we display j/ > only. Clearly 
the unquenched curve describes the data best both in x- 
and y-direction. The onset of Etk is responsible for the 
deviation between data and MM for a; > 10 in Fig. [51 
(left). In Fig.jHlthe scaling of constant ii^Vt is tested for 
these /i-values, where we have again chosen 77^ 3> 100 to 
effectively quench. Both data sets nicely follow the same 
curve, the fits for a agree within errorbars. 

To conclude we have shown that our unquenched 
SU{2) lattice data are quantitatively very well described 
by complex MM predictions over a wide range of chemi- 
cal potentials and masses. We find agreement with e^PT 
by confirming the scaling F^V^^, with a consistent value 
for F^/Yj ~ 0.37 for all our data. This further indicates 



FIG. 6: The scaling of ^4Vi =const.: Fi =6*, /ii =0.2 (red) 
vs. V2 = 8*, =0.1125 (green), both at ma = 20 for y > Q. 

a MM-exPT equivalence for SU{2) at 7^ 0. 
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